%clear all; clc;close all;
%Updated: 11/30/2021

%plot policy functions and value functions
a   =load('agrid2.txt');
b   =load('bgrid2.txt');

vcon0=load('gvcon.txt');
vaut0=load('gvaut.txt');

p0  =load('gp.txt');
h0  =load('gh.txt');
ibp0=load('gibp.txt');
ct0 =load('gct.txt');
ca0 =load('gca.txt');
gn0 =load('ggn.txt');
tau0=load('gtau.txt');
w0  =load('gw.txt');
def0=load('gdef1.txt');
q0  =load('gq.txt');

paut0  =load('gpaut.txt');
haut0  =load('ghaut.txt');
ctaut0 =load('gctaut.txt');
gnaut0 =load('ggnaut.txt');
tauaut0=load('gtauaut.txt');
waut0  =load('gwaut.txt');

na=rows(a);
nk=1;
ns=na/nk;
nb=rows(b);

alpha  = 0.75;
pm     = 1;
gamma  = 0.43;
omegac = 0.3;
mu     = 1;
lambda = 0.184;
r      = 0.02;

atrans = ones(nb,1)*a';

for i=1:na;
    for j=1:nb;
        vcon(i,j)  =vcon0((i-1)*nb + j);
        
        p(i,j)  =p0((i-1)*nb + j);
        h(i,j)  =h0((i-1)*nb + j);
        ibp(i,j)=ibp0((i-1)*nb + j);
        ct(i,j) =ct0((i-1)*nb + j);
        ca(i,j) =ca0((i-1)*nb + j);
        gn(i,j) =gn0((i-1)*nb + j);       
        tau(i,j)=tau0((i-1)*nb + j);       
        w(i,j)  =w0((i-1)*nb + j);   
        def(i,j)=def0((i-1)*nb + j);   
        q(i,j)  =q0((i-1)*nb + j);   
        rec(i,j)=w(i,j)*h(i,j)*tau(i,j);
    end;
end;

bp = b(ibp);

for i=1:na;
    vaut(i)  =vaut0(i);
    
    paut(i)  =paut0(i);
    haut(i)  =haut0(i);
    ctaut(i) =ctaut0(i);
    gnaut(i) =gnaut0(i);       
    tauaut(i)=tauaut0(i);       
    waut(i)  =waut0(i);   
    recaut(i)=waut(i)*haut(i)*tauaut(i);
	wwaut(i) = alpha*paut(i).*haut(i).^(alpha-1); 
end;

cn       = h.^alpha-gn;
cnaut    = haut.^alpha-gnaut;
const    = (omegac.*ct.^(-mu)+(1-omegac).*cn.^(-mu)).^(-1/mu);
constaut = (omegac.*ctaut.^(-mu)+(1-omegac).*cnaut.^(-mu)).^(-1/mu);

for i=1:na;
    temp1 = find(ct(i,:) > 0);
    temp2 = find(cn(i,:) > 0);
    temp3 = find(p(i,:) > 0);

    iminct(i) = temp1(1);
    imincn(i) = temp2(1);
    iminp(i) = temp3(1);
    
    imin(i)  = max(iminct(i),max(imincn(i),iminp(i)));
    
    v(i,1:imin(i)-1)  = NaN; 
    p(i,1:imin(i)-1)  = NaN; 
    bp(i,1:imin(i)-1) = NaN; 
    ibp(i,1:imin(i)-1)= NaN; 
    ctE(i,1:imin(i)-1) = NaN; 
    ctU(i,1:imin(i)-1) = NaN; 
    cnE(i,1:imin(i)-1) = NaN; 
    cnU(i,1:imin(i)-1) = NaN; 
    constE(i,1:imin(i)-1) = NaN;     
    constU(i,1:imin(i)-1) = NaN;     
    gn(i,1:imin(i)-1) = NaN; 
    h(i,1:imin(i)-1)  = NaN; 
    ca(i,1:imin(i)-1)  = NaN; 
    tau(i,1:imin(i)-1)= NaN; 
	rec(i,1:imin(i)-1)= NaN;
	ww(i,1:imin(i)-1)= NaN;
	q(i,1:imin(i)-1)  = NaN;
  
end;

for i=1:ns;
for j=1:nk;        
    vcon2(i,j,:)  = vcon((i-1)*nk+j,:); 
    p2(i,j,:)  = p((i-1)*nk+j,:);  
    bp2(i,j,:) = bp((i-1)*nk+j,:);  
    ibp2(i,j,:)= ibp((i-1)*nk+j,:);  
    ct2(i,j,:) = ct((i-1)*nk+j,:);  
    cn2(i,j,:) = cn((i-1)*nk+j,:);  
    const2(i,j,:) = const((i-1)*nk+j,:);         
    gn2(i,j,:) = gn((i-1)*nk+j,:);  
	ca2(i,j,:) = ca((i-1)*nk+j,:);  
    h2(i,j,:)  = h((i-1)*nk+j,:);   
    tau2(i,j,:)= tau((i-1)*nk+j,:);  
	rec2(i,j,:)= rec((i-1)*nk+j,:); 
    def2(i,j,:)= def((i-1)*nk+j,:);  
	q2(i,j,:)  = q((i-1)*nk+j,:); 
	ww2(i,j,:) = alpha*p2(i,j,:).*h2(i,j,:).^(alpha-1); 
    

    vaut2(i,j)  = vaut((i-1)*nk+j); 
    paut2(i,j)  = paut((i-1)*nk+j);  
    ctaut2(i,j) = ctaut((i-1)*nk+j);  
    cnaut2(i,j) = cnaut((i-1)*nk+j);  
    constaut2(i,j) = constaut((i-1)*nk+j);          
    gnaut2(i,j) = gnaut((i-1)*nk+j);  
    haut2(i,j)  = haut((i-1)*nk+j);  
    tauaut2(i,j)= tauaut((i-1)*nk+j);  
	wwaut2(i,j) = alpha*paut2(i,j).*haut2(i,j).^(alpha-1); 
end;
end;

nettax2    = p2.*gn2-ca2;           %net taxes
nettaxaut2 = paut2.*gnaut2;         %net taxes

fiscal2   = tau2 - p2.*gn2;           %fiscal surplus
fiscalaut2=tauaut2-paut2.*gnaut2;     %fiscal surplus

nfiscal2   = nettax2 - p2.*gn2;           %fiscal surplus net of transfers
nfiscalaut2=nettaxaut2-paut2.*gnaut2;     %fiscal surplus net of transfers

lowa = round(2*ns/6);
meda = round(ns/2);
hgha = round(4*ns/6);

b =-b/(lambda+r);  %change of sign
a = exp(a);

%put together policy fcns of repayment and autarky
mat_ones = ones(size(def));
risk_free_debt=0;

if (risk_free_debt==1);
    p3  =p;
    h3  =h;
    ct3=ct;
    cn3=cn;
    ca3 =ca;
    gn3 =gn;       
    tau3=tau;       
    w3  =w;   
    rec3=rec;
else;
    p3  =def.*p+(mat_ones-def).*(paut'*ones(1,nb));
    h3  =def.*h+(mat_ones-def).*(haut'*ones(1,nb));
    ct3=def.*ct+(mat_ones-def).*(ctaut'*ones(1,nb));
    cn3=def.*cn+(mat_ones-def).*(cnaut'*ones(1,nb));
    ca3 =def.*ca+(mat_ones-def).*0;
    gn3 =def.*gn+(mat_ones-def).*(gnaut'*ones(1,nb));       
    tau3=def.*tau+(mat_ones-def).*(tauaut'*ones(1,nb));       
    w3  =def.*w+(mat_ones-def).*(waut'*ones(1,nb));   
    rec3=def.*rec+(mat_ones-def).*(recaut'*ones(1,nb));
end;

[i1,i2] = find(def>0);
q3  =NaN; bp3 =NaN; sp3 =NaN; gn3 = NaN; h3=NaN; p3=NaN;

for i=1:na;
for j=1:nb;
    if (def(i,j)>0); q3(i,j)=q(i,j);bp3(i,j)=bp(i,j);sp3(i,j)=100*((1+q3(i,j).^(-1)-lambda)./(1+r)-1);p3=p;h3=h;gn3=gn;tau3=tau;w3=w;rec3=rec;
    else
        q3(i,j)=NaN;bp3(i,j)=NaN;sp3(i,j)=NaN;p3=NaN;h3=NaN;gn3=NaN;tau3=NaN;w3=NaN;rec3=NaN;
    end;
end;
end;

nettax3    = p3.*gn3-ca3;           %net taxes
fiscal3    = tau3 - p3.*gn3;        %fiscal surplus
nfiscal3   = nettax3 - p3.*gn3;     %fiscal surplus net of transfers

%value functions: repayment vs autarky
vcon21(:,:) = vcon2(:,1,:);
vcon22(:,:) = vcon2(:,1,:);
vaut21(:,:) = vaut2(:,1);
vaut22(:,:) = vaut2(:,1);

%income index =11 (mean income)

%debt indices
%default IR: 
%=122 mean debt,        value grid = -0.2261
%=196 0.5*mean debt,    value grid = -0.1136
%=159 0.75*mean debt,   value grid = -0.1698

%extra info: bss = 108, value grid = -0.2474

%risk-free debt IR:
%%=1888 mean debt,      value grid = -0.226
%=1944 0.5*mean debt,   value grid = -0.114
%=1916 0.75*mean debt,  value grid = -0.17
%=1774 2*meandebt       value grid = -0.454
%=1877 1.1*mean debt    value grid = -0.248
%=1973 0.25*mean debt,  value grid = -0.056
%=2001 zero,            value grid = 0


figure;
plot(b(1:nb),vcon21(lowa,:),'-r',b(1:nb),vcon21(meda,:),':r',b(1:nb),...
    vaut21(lowa,:)*ones(size(b)),'-b',b(1:nb),vaut21(meda)*ones(size(b)),':b','LineWidth',1.5);
xlabel('b');title('value functions');
ylim([-200 -75]);legend('rep low y^T','rep high y^T','aut low y^T','aut high y^T')

%% Plot variables as function of yT for two debt levels
%optimal policies and prices: repayment vs autarky

if strcmp(label_graph,'policy_fcns');

NameArray = {'LineStyle','Color','LineWidth'};
NameArray0 = {'LineStyle','Color','LineWidth'};
ValueArray0 = {'--','-';[1  0.4  0.4],[0.078  0.1686  0.549];1.5,1.5}';
ValueArray = {'--','-';[1  0.4  0.4],[0.078  0.1686  0.549];1.5,1.5}';
ValueArray2 = {'--','-','--','-';[1  0.4  0.4],[0.078  0.1686  0.549],[1  0.4  0.4],[0.078  0.1686  0.549];1.5,1.5,1.5,1.5}';
ValueArray4 = {'--','--','-','-';'red','blue','red','blue';1.5,1.5,1.5,1.5}';
ValueArray5 = {':',':';[1  0.4  0.4],[0.078  0.1686  0.549];1.5,1.5}';

endb = 2;
lowb = 1916;
hghb = 1888;

nameFig = char([char(label_model) '_' char(label_graph) '_yt']);
figure('name',nameFig);
subplot(5,2,1);P = plot(a,100*squeeze(1-h3(:,lowb)),a,100*squeeze(1-h3(:,hghb)));
title('unemployment','Fontweight','normal');
AX1 = gca;box on;ylabel('%');
set(AX1,'YTick',0:10:20);
axis(AX1, [a(1) a(end) -2 20]);
set(P,NameArray,ValueArray);
subplot(5,2,2);P = plot(a,squeeze(p3(:,lowb)),a,squeeze(p3(:,hghb)));
title('p^{N}','Fontweight','normal'); 
%lgd=legend('low y^T','high y^T');lgd.Location='SouthEast';
legend('low debt','high debt','Location','SouthEast','Orientation','vertical');
AX1 = gca;box on;
set(AX1,'YTick',3.5:0.5:4.5);
axis(AX1, [a(1) a(end) 3.5 4.5]);
set(P,NameArray,ValueArray);
subplot(5,2,3);P = plot(a,squeeze(gn3(:,lowb)),a,squeeze(gn3(:,hghb)));
title('g^{N}','Fontweight','normal'); 
AX1 = gca;box on;
set(AX1,'YTick',0.1:0.1:0.3);
axis(AX1, [a(1) a(end) 0.1 0.3]);
set(P,NameArray,ValueArray);
subplot(5,2,4);P = plot(a,-squeeze(bp3(:,lowb))/(lambda+r),a,-squeeze(bp3(:,hghb))/(lambda+r));
title('debt','Fontweight','normal');
AX1 = gca;box on;
set(AX1,'YTick',0.6:0.3:1.2);
axis(AX1, [a(1) a(end) 0.6 1.2]);
set(P,NameArray0,ValueArray0);
subplot(5,2,5);P = plot(a,squeeze(cn3(:,lowb)),a,squeeze(cn3(:,hghb)));
title('c^{N}','Fontweight','normal');
AX1 = gca;box on;
set(AX1,'YTick',0.5:0.2:0.9);
axis(AX1, [a(1) a(end) 0.5 0.9]);
set(P,NameArray,ValueArray);
subplot(5,2,6);P = plot(a,squeeze(ct3(:,lowb)),a,squeeze(ct3(:,hghb)));
title('c^{T}','Fontweight','normal');
AX1 = gca;box on;
set(AX1,'YTick',0.7:0.2:1.1);
axis(AX1, [a(1) a(end) 0.7 1.1]);
set(P,NameArray,ValueArray);
subplot(5,2,7);P = plot(a,squeeze(w3(:,lowb)),a,squeeze(w3(:,hghb)));
title('wages','Fontweight','normal');
AX1 = gca;box on;
set(AX1,'YTick',2.5:0.5:3.5);
axis(AX1, [a(1) a(end) 2.5 3.5]);
set(P,NameArray,ValueArray);
%version 1: with bond spreads
subplot(5,2,8);P = plot(a,squeeze(sp3(:,lowb)),a,squeeze(sp3(:,hghb)));
title('spreads','Fontweight','normal');
set(P,NameArray0,ValueArray0);ylabel('%');
AX1 = gca;box on;
set(AX1,'YTick',0:4:8);
axis(AX1, [a(1) a(end) -0.5 8]);
%version 2: with current account
subplot(5,2,9);P = plot(a,-squeeze(ca3(:,lowb)),a,-squeeze(ca3(:,hghb)));
title('CA balance','Fontweight','normal');
AX1 = gca;box on;
set(AX1,'YTick',-0.1:0.1:0.2);
axis(AX1, [a(1) a(end) -0.1 0.2]);
set(P,NameArray0,ValueArray0);
%version 3: with (net) taxes
subplot(5,2,10);P = plot(a,squeeze(tau3(:,lowb)),a,squeeze(tau3(:,hghb)));
set(P,NameArray,ValueArray);
hold on;
PP = plot(a,squeeze(nettax3(:,lowb)),a,squeeze(nettax3(:,hghb)));
title('(net) taxes','Fontweight','normal');
AX1 = gca;box on;
set(AX1,'YTick',0.6:0.3:1.2);
axis(AX1, [a(1) a(end) 0.6 1.2]);
set(PP,NameArray,ValueArray5);
saveas(gcf, fullfile(plot_path, nameFig),format_chart);

%% Plot variables as function of debt for two yT levels
hr2  = h2;
cnr2 = cn2;
gnr2 = gn2;
ctr2 = ct2;
pr2  = p2;
taur2= tau2;
nettaxr2= nettax2;
fiscalr2  = fiscal2;
nfiscalr2 = nfiscal2;
bpr2 = bp2;
wwr2 = ww2;
qr2  = q2;
spr2 = 100*((1+qr2.^(-1)-lambda)./(1+r)-1);
car2 = ca2;

%optimal policies and prices: repayment vs autarky
NameArray = {'LineStyle','Color','LineWidth'};
ValueArray = {'--','--','-','-';[1  0.4  0.4],[1  0.4  0.4],[0.078  0.1686  0.549],[0.078  0.1686  0.549];1.5,1.5,1.5,1.5}';
ValueArray2 = {'--','-','--','-';[1  0.4  0.4],[0.078  0.1686  0.549],[1  0.4  0.4],[0.078  0.1686  0.549];1.5,1.5,1.5,1.5}';
ValueArray4 = {'--','--','-','-';'red','blue','red','blue';1.5,1.5,1.5,1.5}';
ValueArray5 = {':',':',':',':';[1  0.4  0.4],[1  0.4  0.4],[0.078  0.1686  0.549],[0.078  0.1686  0.549];1.5,1.5,1.5,1.5}';

NameArray0 = {'LineStyle','Color','LineWidth'};
ValueArray0 = {'--','-';[1  0.4  0.4],[0.078  0.1686  0.549];1.5,1.5}';

[ilowa_aut]=find(def2(lowa,1,:)==0);
[ilowa_rep]=find(def2(lowa,1,:)==1);
hr2(lowa,1,ilowa_aut) = NaN;
ha2(lowa,1,1:nb) = haut2(lowa,1);
ha2(lowa,1,ilowa_rep) = NaN;

[ihgha_aut]=find(def2(hgha,1,:)==0);
[ihgha_rep]=find(def2(hgha,1,:)==1);
hr2(hgha,1,ihgha_aut) = NaN;
ha2(hgha,1,1:nb) = haut2(hgha,1);
ha2(hgha,1,ihgha_rep) = NaN;

pr2(lowa,1,ilowa_aut) = NaN;
pa2(lowa,1,1:nb) = paut2(lowa,1);
pa2(lowa,1,ilowa_rep) = NaN;

pr2(hgha,1,ihgha_aut) = NaN;
pa2(hgha,1,1:nb) = paut2(hgha,1);
pa2(hgha,1,ihgha_rep) = NaN;

gnr2(lowa,1,ilowa_aut) = NaN;
gna2(lowa,1,1:nb) = gnaut2(lowa,1);
gna2(lowa,1,ilowa_rep) = NaN;

gnr2(hgha,1,ihgha_aut) = NaN;
gna2(hgha,1,1:nb) = gnaut2(hgha,1);
gna2(hgha,1,ihgha_rep) = NaN;

bpr2(lowa,1,ilowa_aut) = NaN;
bpr2(hgha,1,ihgha_aut) = NaN;

cnr2(lowa,1,ilowa_aut) = NaN;
cna2(lowa,1,1:nb) = cnaut2(lowa,1);
cna2(lowa,1,ilowa_rep) = NaN;

cnr2(hgha,1,ihgha_aut) = NaN;
cna2(hgha,1,1:nb) = cnaut2(hgha,1);
cna2(hgha,1,ihgha_rep) = NaN;

ctr2(lowa,1,ilowa_aut) = NaN;
cta2(lowa,1,1:nb) = ctaut2(lowa,1);
cta2(lowa,1,ilowa_rep) = NaN;

ctr2(hgha,1,ihgha_aut) = NaN;
cta2(hgha,1,1:nb) = ctaut2(hgha,1);
cta2(hgha,1,ihgha_rep) = NaN;

wwr2(lowa,1,ilowa_aut) = NaN;
wwa2(lowa,1,1:nb) = wwaut2(lowa,1);
wwa2(lowa,1,ilowa_rep) = NaN;

wwr2(hgha,1,ihgha_aut) = NaN;
wwa2(hgha,1,1:nb) = wwaut2(hgha,1);
wwa2(hgha,1,ihgha_rep) = NaN;

taur2(lowa,1,ilowa_aut) = NaN;
taua2(lowa,1,1:nb) = tauaut2(lowa,1);
taua2(lowa,1,ilowa_rep) = NaN;

taur2(hgha,1,ihgha_aut) = NaN;
taua2(hgha,1,1:nb) = tauaut2(hgha,1);
taua2(hgha,1,ihgha_rep) = NaN;

nettaxr2(lowa,1,ilowa_aut) = NaN;
nettaxa2(lowa,1,1:nb) = nettaxaut2(lowa,1);
nettaxa2(lowa,1,ilowa_rep) = NaN;

nettaxr2(hgha,1,ihgha_aut) = NaN;
nettaxa2(hgha,1,1:nb) = nettaxaut2(hgha,1);
nettaxa2(hgha,1,ihgha_rep) = NaN;

fiscalr2(lowa,1,ilowa_aut) = NaN;
fiscala2(lowa,1,1:nb) = fiscalaut2(lowa,1);
fiscala2(lowa,1,ilowa_rep) = NaN;

fiscalr2(hgha,1,ihgha_aut) = NaN;
fiscala2(hgha,1,1:nb) = fiscalaut2(hgha,1);
fiscala2(hgha,1,ihgha_rep) = NaN;

nfiscalr2(lowa,1,ilowa_aut) = NaN;
nfiscala2(lowa,1,1:nb) = nfiscalaut2(lowa,1);
nfiscala2(lowa,1,ilowa_rep) = NaN;

nfiscalr2(hgha,1,ihgha_aut) = NaN;
nfiscala2(hgha,1,1:nb) = nfiscalaut2(hgha,1);
nfiscala2(hgha,1,ihgha_rep) = NaN;

car2(lowa,1,ilowa_aut) = NaN;
car2(hgha,1,ihgha_aut) = NaN;

qr2(lowa,1,ilowa_aut) = NaN;
qr2(hgha,1,ihgha_aut) = NaN;

spr2(lowa,1,ilowa_aut) = NaN;
spr2(hgha,1,ihgha_aut) = NaN;

endb = 2;

nameFig = char([char(label_model) '_' char(label_graph) '_debt']);
figure('name',nameFig);
subplot(5,2,1);P = plot(b,100*squeeze(1-hr2(lowa,1,:)),b,100*squeeze(1-ha2(lowa,1,:)),b,100*squeeze(1-hr2(hgha,1,:)),b,100*squeeze(1-ha2(hgha,1,:)));
title('unemployment','Fontweight','normal');
AX1 = gca;box on;ylabel('%');
set(AX1,'YTick',0:10:20);
axis(AX1, [min(b) endb -1 20]);
set(P,NameArray,ValueArray);
subplot(5,2,2);P = plot(b,squeeze(pr2(lowa,1,:)),b,squeeze(pr2(hgha,1,:)),b,squeeze(pa2(lowa,1,:)),b,squeeze(pa2(hgha,1,:)));
title('p^{N}','Fontweight','normal'); 
legend('low y^T','high y^T','Location','NorthEast','Orientation','vertical');
AX1 = gca;box on;
set(AX1,'YTick',3.5:2:7.5);
axis(AX1, [min(b) endb 3.5 7.5]);
set(P,NameArray,ValueArray2);
subplot(5,2,3);P = plot(b,squeeze(gnr2(lowa,1,:)),b,squeeze(gna2(lowa,1,:)),b,squeeze(gnr2(hgha,1,:)),b,squeeze(gna2(hgha,1,:)));
title('g^{N}','Fontweight','normal'); 
AX1 = gca;box on;
set(AX1,'YTick',0.18:0.05:0.28);
axis(AX1, [min(b) endb 0.18 0.28]);
set(P,NameArray,ValueArray);
subplot(5,2,4);P = plot(b,-squeeze(bpr2(lowa,1,:))/(lambda+r),b,-squeeze(bpr2(hgha,1,:))/(lambda+r));
title('debt','Fontweight','normal');
AX1 = gca;box on;
set(AX1,'YTick',-14:8:2);
axis(AX1, [min(b) endb -14 2]);
set(P,NameArray0,ValueArray0);
subplot(5,2,5);P = plot(b,squeeze(cnr2(lowa,1,:)),b,squeeze(cna2(lowa,1,:)),b,squeeze(cnr2(hgha,1,:)),b,squeeze(cna2(hgha,1,:)));
title('c^{N}','Fontweight','normal');
AX1 = gca;box on;
set(AX1,'YTick',0.6:0.1:0.9);
axis(AX1, [min(b) endb 0.6 0.9]);
set(P,NameArray,ValueArray);
subplot(5,2,6);P = plot(b,squeeze(ctr2(lowa,1,:)),b,squeeze(cta2(lowa,1,:)),b,squeeze(ctr2(hgha,1,:)),b,squeeze(cta2(hgha,1,:)));
title('c^{T}','Fontweight','normal');
AX1 = gca;box on;
set(AX1,'YTick',0.8:0.2:1.4);
axis(AX1, [min(b) endb 0.8 1.4]);
set(P,NameArray,ValueArray);
subplot(5,2,7);P = plot(b,squeeze(wwr2(lowa,1,:)),b,squeeze(wwa2(lowa,1,:)),b,squeeze(wwr2(hgha,1,:)),b,squeeze(wwa2(hgha,1,:)));
title('wages','Fontweight','normal');xlim([min(b) max(b)]);
AX1 = gca;box on;
set(AX1,'YTick',3:1:6);
axis(AX1, [min(b) endb 3 6]);
set(P,NameArray,ValueArray);
%version 1: with bond spreads
subplot(5,2,8);P = plot(b,squeeze(spr2(lowa,1,:)),b,squeeze(spr2(hgha,1,:)));
title('spreads','Fontweight','normal');
set(P,NameArray0,ValueArray0);ylabel('%');
AX1 = gca;box on;
set(AX1,'YTick',0:20:40);
axis(AX1, [min(b) endb -2 40]);
%version 2: with current account
subplot(5,2,9);P = plot(b,-squeeze(car2(lowa,1,:)),b,-squeeze(car2(hgha,1,:)));
title('CA balance','Fontweight','normal');
AX1 = gca;box on;
set(AX1,'YTick',-0.3:0.3:0.3);
axis(AX1, [min(b) endb -0.3 0.3]);
set(P,NameArray0,ValueArray0);
%version 3: with (net) taxes
subplot(5,2,10);P = plot(b,squeeze(taur2(lowa,1,:)),b,squeeze(taua2(lowa,1,:)),b,squeeze(taur2(hgha,1,:)),b,squeeze(taua2(hgha,1,:)));
set(P,NameArray,ValueArray);
hold on;
PP = plot(b,squeeze(nettaxr2(lowa,1,:)),b,squeeze(nettaxa2(lowa,1,:)),b,squeeze(nettaxr2(hgha,1,:)),b,squeeze(nettaxa2(hgha,1,:)));
title('(net) taxes','Fontweight','normal');
AX1 = gca;box on;
set(AX1,'YTick',0.25:0.5:1.75);
axis(AX1, [min(b) endb 0.25 1.75]);
set(PP,NameArray,ValueArray5);
saveas(gcf, fullfile(plot_path, nameFig),format_chart);
end;
